
lambda = 1;
N = 8;
d = lambda / 2;

mu = 0.05;
ML = 1000;

h = modem.pskmod('M', 4);              % Modulator object
msg = randi([0 3],ML,1);               % Modulating message
x = (modulate(h,msg) * ones(1,N)).';
x = x ./ N;

% FROM HERE, IMPLEMENTATION FROM KOEN
% M  = zeros(N,N);                    % Generate matrix of form:
%                                     %   0           ...   n_0-n_{N-1};   
% for r=1:(N)                         %   .           0                ; 
%    M(r,:)=(r-1:(-1):-N+r);          %  n_{N-1}-n_0  ...             0;]
% end
% 
% nv           = (0:1:N-1).';
% % Define phi and phi_prime (as in paper)
% phi         = @(theta) (1j*2*pi*d*sind(theta))/lambda;
% phi_p       = @(theta) (1j*2*pi*d*cosd(theta))/lambda;
% 
% R2 = 1;
% a_prime         = @(theta_k,xv) real( (xv'*( phi_p(theta_k).*M.*exp(M.*phi(theta_k)) )*xv) );  
% theta_k_next    = @(theta_k,xv) theta_k - mu*(abs( exp(phi(theta_k)*nv)'*xv )^2 - R2) * a_prime(theta_k,xv);
% END HERE

k=1;
theta = -14;
values = zeros(ML - 1, 1);
while k < ML,
   theta_old = theta;
   w = calc_w(theta, lambda, N);
   y = (w' * x(:, k));
   BB = calc_bb(theta, lambda, N);
   
   values(k) = theta;
   theta = real(theta_old - mu*(abs(y)^2-1)*(x(:, k)'*BB*x(:, k)))
   %theta = theta_old - mu*(abs(y)^2-1)*a_prime(theta_old, x(:, k))
   %theta = theta_k_next(theta_old, x(:, k));
   k=k+1;
   
end

plot(values);